Defect Motion and Lattice Pinning Barrier in Josephson-Junction Ladders 
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We study motion of domain wall defects in a fully frustrated Josephson-junction ladder system, 
driven by small applied currents. For small system sizes, the energy barrier Eb to the defect motion 
is computed analytically via symmetry and topological considerations. More generally, we perform 
numerical simulations directly on the equations of motion, based on the resistively-shunted junction 
model, to study the dynamics of defects, varying the system size. Coherent motion of domain walls 
is observed for large system sizes. In the thermodynamical limit, we find Eb = 0.1827 in units of 
the Josephson coupling energy. 

PACS numbers: 74.50.+r, 03.75.Lm, 74.81. Fa, 74.25.Qt 



I. INTRODUCTION 

Two-dimensional (2D) arrays of Josephson junctions 
are of interest in various fields of fundamental classi- 
cal and quantum physics. In the simplest case, they 
provide an experimental realization of the XY model; 
in particular, applying a magnetic field introduces frus- 
tration, measured by the flux per plaquette in units 
of the flux quantumi The corresponding vortices in- 
duced by the field tend to form a regular flux lattice, 
thus lowering the free energy, and result in interplay 
with the underlying lattice periodicity. This gives rise 
to commensurate-incommensurate effects and leads to a 
rich variety of physics, including first-order and double 
transitions, reentrance, glassy behavior, quantum transi- 
tions, topological quantization, dynamic transitions and 
resonance, etc^£ In these phenomena vortex configura- 
tions and dynamics play crucial roles, driving transitions 
and governing transport properties. Here one interesting 
question arises when an extra vortex is added into the 
system. While the vortex in general sits on a plaquette 
with minimum energy, which is separated by the poten- 
tial barrier set by the underlying lattice structure, it may 
be driven into motion by applying currents, as it is ex- 
erted by the "Lorentz force" in the transverse direction, 
and accordingly generates non- vanishing voltage. Indeed, 
the voltage measurement in recent dynamic simulations, 4 
performed in the presence of external currents, has given 
the pinning energy barriers as well as the critical cur- 
rents, which agree fully with experimental results^ thus 
resolved the long-standing discrepancy in the frustrated 
case. 

This paper focuses on the vortex dynamics in ladders 
of Josephson junctions, which provides the simplest sys- 
tem for probing the frustration effects: Those studied in 
existing literature include the vortex configuration and 
the critical current, depending on the frustration the 
vortex-vortex interaction decaying exponentially^ quan- 
tum effects^, and resonance^ Note the vast difference 
from the 2D system, especially, in the vortex interaction, 



which is expected to affect significantly the dynamics of 
a vortex in a background of other vortices, i.e., in a frus- 
trated system. In particular domain walls in a ladder 
system assume the simple form of point defects, the dy- 
namics of which is convenient to probe. We thus consider 
the domain wall defects created by adding an extra vor- 
tex in a fully frustrated ladder and examine their motion 
driven by external currents. In small systems, the sym- 
metry argument and topological constraints allow one to 
computed analytically the energy barrier. More gener- 
ally, the defect motion, driven by uniform currents, is in- 
vestigated by means of dynamical simulations performed 
directly on the equations of motion. The resulting value 
of the energy barrier is found consistent with the ana- 
lytical one obtained for small systems. Also observed is 
the defect motion, either sequential or simultaneous, de- 
pending on the size and the initial configuration. Such 
characteristics are attributed to the distance-dependent 
interaction between defects and the underlying lattice ge- 
ometry. 

There are five sections in this paper: Section II intro- 
duces the model system whereas Sec. Ill is devoted to 
the analytical calculations of the energy barrier to the 
defect motion in small systems. In Sec. IV, we describe 
the numerical simulations performed on the equations of 
motion in the presence of uniform driving currents, and 
present the results. The current-voltage (IV) character- 
istics and the energy function are computed, which in 
turn give the critical current and the pinning energy bar- 
rier for various system sizes. Finally, a brief summary is 
given in Sec. V. 



II. MODEL SYSTEM 

We consider a ladder of Josephson junctions made of 
2L superconducting grains weakly coupled to their near- 
est neighbors, the schematic diagram of which is shown 
in Fig.^ The grains are located at sites i = (x, y), where 
x runs from 1 to L (in the leg direction) and the label 
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FIG. 1: Schematic notation for a Josephson-j unction ladder. 
Each superconducting grain, denoted by a square, is charac- 
terized by the phase 4> x , y of the superconducting order param- 
eter. The symbol ® denotes the flux per plaquette from an 
external transverse magnetic field. The extra plaquette on the 
right hand side represents the periodic boundary conditions. 



y (— 1,2) describes respectively the lower and upper legs 
of the ladder. Each grain is characterized by the local 
condensed wave function or the order parameter: 



*i = |*i|e 



icf>i 



(i) 



where the local superconducting fluid density is as- 
sumed to be constant at low temperatures. Accordingly, 
relevant fluctuations come from the phases (\>i and the 
Hamiltonian of the system in the presence of the external 
field is simply given by the sum of the nearest neighbor- 
ing pair energies 



H = -E 



E 



cos 



A; 



(2) 



where Ej is the coupling constant between the grains, 
(i, j) represents nearest neighboring pairs, and the bond 
angle A+j is given by the line integral of the vector po- 
tential: 



A; 



2vr 



A-dl 



(3) 



with the flux quantum $o = nhc/e. In the Landau gauge, 
the components of the vector potential A(x, y) are given 
by 



A x (x,y)=0 and A y (x, y) = $x, 



(4) 



where $ is the magnetic flux per plaquette and x is the 
position along the leg direction. 

For the ladder in Fig. ^ the Hamiltonian reduces to 



H = - Ej^2 cos ( 



Ej y ,cos^ x ,i 



•x,2 - 2?r/x) 



(5) 



where / = $/$o measures the frustration of the system. 
In the fully frustrated case (/ = 1/2), which is our main 
concern in this work, one every other site is occupied by 
a single vortex. 

We now add or remove one vortex; this creates topolog- 
ical defects (domain walls) that affect the ground state. 
A typical vortex configuration in this case is displayed 



FIG. 2: Vortex configuration in the presence of an extra vor- 
tex in the fully frustrated ladder of size L = 16. Filled circles 
represent vortices. 



in Fig. |21 The extra vortex can move through the peri- 
odic potential produced by the lattice structure when it is 
subject to a perpendicular current. An estimation of the 
corresponding lattice pinning barrier is then made each 
time this extra vortex crosses the barrier. Note that the 
periodic potential is in general modulated significantly 
by other (underlying) vortices present in the system with 
/ = 1/2, resulting in the barrier strikingly different from 
that in the unfrustrated system (/ = 0). 



III. ANALYTIC CALCULATIONS 

For convenience, we choose new gauge invariant phases 
that simplify the Hamiltonian and the current distribu- 
tion in the system. Let 9 X and ip x denote the following 
phase differences between the grains: 



i> x +i,2 — 4>x,2 



: fj. 



»x,2 



2tt/x. 



(0) 



It is easy to see, by symmetry and energy considerations, 
that the phase differences 9 X and 6' x are opposite to each 
other. 6 Indeed, the sum of these phases around each pla- 
quette is constrained topologically by the flux or frustra- 
tion / and the (integer) vortex number n x : 



<Px+i +<Px = 2Tr(n x - /) = 2Tvq x 



(7) 



where q x is the (fractional) vortex charge, and the Hamil- 
tonian simply reads 

H = —Ej [cos 9 X + cos 6' x + cos ip x ] 



2 cos 



<Px+l - fx 



+ n(n x - /) 



COS (fx 



(8) 



Then the condition 9' x — —9 X decouples the phases be- 
tween the transverse directions and leads to a solution 
that minimizes this Hamiltonian. Using the current con- 
servation laws, we can write a set of L equations for 9 X 
and tfx at every node of the lattice: 



sin = sm0 x+ i - sin</j x+ i 
with the boundary conditions 

9 x +l = 9 X and ip x+L = ip x . 



(9) 
(10) 
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FIG. 3: Configuration of a two-piaquette system, with the 
phase difference iabeied aiong each fink. Filied circies repre- 
sent vortices. 
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FIG. 4: Phase configuration of a Josephson-j unction ladder 
of L = 8 plaquettes. Periodic boundary conditions are em- 
ployed. 



The barrier energy Eb for a vortex moving along a lad- 
der can be computed exactly on one simple example. In 
Fig. |3| we consider two plaquettes under closed bound- 
ary conditions and a single vortex in the system. In the 
notation of Fig. |3 the equations for the phases (0, 9' , p) 
in case (a) are given by 



30- 
W 



sm v = sm f 



- smip, 



(11) 



which yields 9 = —9' = (tt + ip)/3. As a function of p, 
the energy 



E(p) = —6 cos 



TT + p 



cos p 



(12) 



has an absolute minimum for p = -—tt/2, which in turn 
leads to = tt/6 and E — — 3y3 « —5.196, and a maxi- 
mum for <p = — tt together with 9 = and E = — 5. On 
the other hand, in case (b), we have 



30- p 
30' + p 



sm f 



— 7T 

n 

sin0 + ship, 



(13) 



the solutions of which are p = tt/2 and 9 = —tt/6 for the 
ground state (with E = — 3\/3) and p = n and 9 = 
for the excited state (with E = —5). The excited states 
in both cases are equivalent since p = —tt = tt (mod 2ir). 
This corresponds to the situation that the system evolves 
from configuration (a) to (b), namely, the instant when 
the vortex is exactly on the rung between the two pla- 
quettes. Accordingly, the energy barrier is simply given 
by 



E B = 3\/3-5 w 0.196. 



(14) 



In the general case, the value p = tt (or — tt) in the 
excited state does not depend on the frustration pa- 
rameter / since it is always a solution of the equation 
dE{p)/dp = with 



hyp) = — A cos I | — 3 cos 



cos p. 



2tt(1 -/) + ¥> 
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-16.3229 



-16.1368 



(15) 



FIG. 5: Three vortex configurations (I, G,M), showing the 
presence of vortices (denoted by filled circles) and domain 
wall defects (denoted by crosses). Also shown is the esti- 
mated energy (in units of Ej) of each configuration. Arrows 
represent the time evolution of the configuration, which has 
been verified by extensive numerical simulations. 



In the following, we accept that p = u corresponds to 
the solution of the excited state in which the vortex is on 
the rung for any given L; this will be checked numerically 
(see Fig. below). 

In a more complicate case, we consider a system with 
L = 8 plaquettes. Figure0]shows the phase configuration 
of an eight-plaquette system at / = 1/2 under periodic 
boundary conditions, in the presence of an extra vortex. 
We consider three possible configurations (I, G, M) shown 
in Fig. [SI where filled circles and crosses represent vor- 
tices and defects (domain walls), respectively. Starting 
from the initial state I and driven by the injected cur- 
rent along the y direction, the system evolves eventually 
to configuration G via a number of intermediate config- 
urations. It subsequently evolves to M and back to G. 
Configurations G and M correspond to the lowest-energy 
state and the high-energy (excited) state, respectively, 
and this evolution pattern repeats with time, which has 
been verified by extensive numerical simulations. 

In the initial configuration I, the stationary phase re- 
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(16) 



In units of the Josephson coupling energy Ej, the energy 
is estimated to be E(I) — —16.3083. For configuration 
G, the phase relations read 



201 - (f>2 + 9l 
20 3 - (f)4 + ip 3 

205 - ^6 + 
20 7 - (f) s + tp 7 



"71", 
-7T, 
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208 



93 
95 
tp 7 

<Pl 



92 
96 



7T 
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7T 



(17) 



which yields the energy E(G) = -16.3229. 

Configuration M describes an intermediate state via 
which the system goes from the state with the occupation 
number (ji 5 = 0, Uq = 1) to that with (715 = 1, rig = 0), 
namely, the vortex moves to the left by one plaquette, 
similarly to the evolution from (b) to (a) in Fig. [21 In 
this case, the vortex numbers in both cells are not well 
defined, but the vortex is said to be "spread" between 
the two plaquettes. To apply Eq. J§J, we further take 
the two plaquettes on both sides of the rung as one unit 
cell. Since the net vortex charge enclosed in this cell (con- 
sisting of the two plaquettes) is zero, the sum of phase 
differences around it also vanishes. We thus have the 
condition 205 + tp 5 + 206 — ip 7 = 27r(n5 — / + tiq — /) = 0, 
with n 5 + n e — 1. The remaining relations are given by 



20i - tp 2 

20 3 - tp 4 

20 7 - p 8 
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(18) 



As addressed already, the vortex sits on the rung in this 
configuration and the phases take the radial direction 
around the center of the rung, thus leading to the phase 
difference ipe — ir along the rung. This is manifested by 
the time evolution of ipe, as shown in the next section 
(see Fig. EJ). We thus set ipe = ir and obtain the energy 
of the configuration: E(M) = -16.1368. 

Together with the result of E(G), we estimate the pin- 
ning barrier according to 



E B ee E(M) - E(G) = 0.1861. 



(19) 



Note that this value, obtained for L = 8, is lower than the 
value 0.19615 in the two-plaquette case (L = 2). We thus 
expect that the energy barrier Eb in the thermodynamic 
limit (L — > 00) has a value still lower than 0.1861. 



IV. NUMERICAL SIMULATIONS 

To evaluate the precise value of the energy barrier for 
various system sizes, we have performed extensive dy- 
namic simulations on the resistively shunted junction 
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FIG. 6: Time evolution of the phase difference across the 
rung (solid line, left vertical axis), together with that of the 
vortex charge q = n — / on the plaquette just right of the 
rung (dotted line, right vertical axis). Time t is given in units 
of h/2ei c R. 



(RSJ) model. The dynamics of the RSJ model, with 
single-junction critical current i c and shunt resistance R, 
is governed by the set of equations of motion for the phase 



2eR dt 



(20) 



where 



A; 



is the gauge-invariant phase 
difference across the junction (ij), and the primed sum- 
mation runs over the nearest neighbors of grain i. The 
system is driven by the current U — l x , y — I{5y,2 — $y,i) 
(applied to grain i), namely, uniform current / is injected 
to and extracted from each grain on the upper (y = 2) 
and lower (y = 1) legs, respectively. Using a modified 
Euler method, we have integrated Eq. I|20|) with the time 
step of size At = 0.05 (in units of h/2ei c R) for a variety 
of ladders up to the system size L = 512. In addition 
to the periodic boundary conditions imposed along the x 
direction, we introduce a 27r phase slip across the whole 
system: 



9L + 1,2 = 01,2 + 7T 
t>L + l.l = 4>1,1 - 7T, 



(21) 



which generates a single extra vortex. 

We first examine how the rung phase difference (p§ 
varies in the vortex motion and plot in Fig. its time 
evolution in the system of eight plaquettes. Also plotted 
is the evolution of the vortex charge q = n — / (with 
n being the vortex number) on the plaquette just right 
of the rung. It is observed that q (or n) changes rather 
abruptly from 1/2 to —1/2 (or from 1 to 0), describing 
the motion of a vortex to the left. In particular, at the 
moment of the change, i.e., when the vortex is located on 
the rung, the phase difference ipe indeed has the value 7r, 
as expected. 
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FIG. 7: Pattern of defect motion in a fully frustrated ladder 
(/ = 1/2), with filled circles and crosses denoting vortices and 
domain walls, respectively, as time goes by (in the direction of 
the arrows). Currents are applied uniformly along the rungs. 
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FIG. 8: IV characteristics for the system size L = 8, 16, and 
24, respectively, from bottom to top. The inset shows the 
critical current I c as a function of L. Current I and voltage 
V are expressed in units of i c and i c R, respectively. 



Figure shows typical motion of defects under the 
driving currents. At first, two defects (i.e., two domain 
walls separating the three neighboring vortices) are next 
to each other, as shown in the first configuration (from 
top to bottom). The distance between the domain wall 
defects grows with time until this distance eventually be- 
comes half the system size (see the second configuration) . 
Then, the defect on the right moves first (changing the 
configuration to the third one) , subsequently followed by 
the motion of the one on the left (resulting in the fourth 
configuration). In the case that there are only a few 
plaquettes (L < 40), this behavior is always observed, 
regardless of the initial distance between the two defects. 
On the other hand, in a system of larger size, two types 
of behavior are observed, depending on initial conditions: 
When the two defects are initially located at nearby sites, 
they move simultaneously through transient states and 
the distance between them does not grow beyond 20 pla- 
quettes. In contrast, two defects distant by more than 20 
plaquettes tend to move sequentially for appropriate ini- 
tial phase configurations. We presume that such size de- 
pendence has its origin in the interaction between defects 



and the underlying periodic lattice geometry. Namely, 
the interaction between two domain wall defects becomes 
vanishingly small as the distance is increased beyond 20 
plaquettes, which may reflect the exponentially decaying 
interaction between vortices^ In this manner the charac- 
teristic interaction between domain walls in a background 
of vortices appears to be exposed. 

In order to estimate the pinning barrier, we compute 
the IV characteristics and the critical current I c , and 
probe their behaviors with the system size L. The voltage 
across the system is given by the ac Josephson relation-ii 
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FIG. 9: Time evolution of energy E(t) in systems of size 
L = 8 (a) and 64 [(b) and (c)]. The uniform driving cur- 
rent I = I C (L) +0.0001 has been applied along each rung, (a) 
and (b) describe the sequential motion while (c) corresponds 
to the simultaneous motion (see the text). For convenience, 
E(t), given in units of E.j, has been shifted such that E = 
corresponds to the minimum. The inset in (c) shows a tran- 
sient behavior: The two peaks merge eventually into one peak 
shown in the main plate. 
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FIG. 10: Pinning energy barrier Eb (in units of Ej) as a 
function of the system size L. Energy barriers corresponding 
to (a) sequential and (b) simultaneous motion of defects are 
displayed. 



respectively. As the defect moves across one plaquette, 
E(t) goes through a maximum corresponding to the ex- 
cited state discussed in Sec. III. The lowest-energy state 
corresponds to configuration G and the maximum one 
to M shown in Fig. [5J As pointed out, the defects can 
move simultaneously for appropriate initial conditions. 
Such simultaneous motion is indeed observed in Fig.|Sfc), 
which reveals the doubling of both the amplitude and the 
period of E(t) (i.e., the energy barrier and the period of 
the defect motion). The two transient states seen in the 
inset of Fig.|^c) indicate that the system possessing two 
defects is not completely coherent in the first stage of the 
dynamics. 

The pinning energy barrier Eb, defined to be the differ- 
ence between the maximum energy E(M) and the min- 
imum one E(G), is thus computed as the system size is 
varied. The size dependence of Eb is then examined and 
shown in Fig. I1UI for (a) sequential and (b) simultaneous 
motion of defects. In the former case, the energy barrier 
is observed to approach the value 



E B = 0.1827 



(23) 



in the thermodynamic limit. This value is slightly be- 
low the one found analytically in the eight-plaquette sys- 
tem, as expected. In the case of simultaneous motion, 
Fig. HUf bl shows that the energy barrier becomes double 
for the system size L > 40. 



V. SUMMARY 



and the resulting characteristics are displayed in Fig. |S| 
for system size L = 8, 16, and 24. Systems larger than 
L = 24 turn out to exhibit the same IV characteristics 
as the case L = 24 and are thus not shown here. It 
is observed that non-zero voltage develops as the driv- 
ing current I is increased beyond a certain value. The 
size-dependence of the corresponding critical current I c 
is plotted in the inset of Fig. [H] which demonstrates that 
I c first reduces with the system size L and saturates to 
a nearly constant value beyond L = 24. In the thermo- 
dynamic limit, I c is shown to approach the value 0.089 
(in units of the single-junction critical current i c ); this is 
close to the value I c « 0.1 at f — 1/2, extracted roughly 
from Fig. 1(a) of Ref.Q. 

In Fig. we display the typical time dependence of the 
energy E(t). With the driving current I = I C (L) + 0.0001 
just above the critical value, the energy is calculated 
through the use of Eq. JSJ). Note in Fig. [5] that (a) and 
(b) correspond to sequential motion of defects for the sys- 
tem size L = 8 (smaller than 40) and 64 (larger than 40), 



We have studied the dynamics of domain wall defects 
created by adding an extra vortex in a fully frustrated 
Josephson-junction ladder. The defects are in general 
pinned by the energy barrier generated by the underlying 
lattice structure and other vortices induced by an exter- 
nal magnetic field or frustration. Making use of the sym- 
metry and topological constraints, we have computed the 
energy barrier Eb in systems of size L < 8. The defects 
may be put to motion by applying currents larger than 
the critical current. The corresponding motion in the 
system, driven by uniform currents just above the criti- 
cal value, has been investigated by means of dynamical 
simulations performed directly on the equations of mo- 
tion. The resulting numerical estimation of Eb — 0.1827 
(in units of the Josephson coupling energy) is fully con- 
sistent with the analytical value obtained from resolu- 
tion of the phases in the eight-plaquette (L = 8) sys- 
tem. In the dynamical study of the system, we have also 
observed that the defects move sequentially in small sys- 
tems (L < 40). On the other hand, in larger systems, the 
domain walls may also display coherent motion, namely, 
they can move simultaneously as well as sequentially, de- 
pending on the initial configurations. Such difference has 
been attributed to the distance-dependent interaction be- 
tween defects and the underlying lattice geometry. 
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